Analyse du séisme de Hebgen Lake (1959)

Conséquences socio-économique de ce séisme

Ce séisme a eu lieu en 1959, dans l'état du Wyoming aux Etats-Unis. L'épicentre se situait en plein milieu du parc de Yellowstone. C'était un séisme de magnitude 7.3 mw, et la faille s'étendait jusqu'à 5km de profondeur.

Conséquences du séisme :

• Ce séisme a eu lieu en pleine forêt, dans un parc naturel, quasiment aucun batiment n'a été endommagé ni détruit, il n'y a donc pas de dommage aux bâtiments.
• Très peu de phénomènes de liquéfaction ont été reportés, ils n'ont donc eu aucun impact significatif sur le paysage.
• Les mouvements gravitaires des ondes sismiques ont cependant été très forts : glissements de terrains, chutes de blocs, coulées de boues...
• L'impact le plus important est aussi l'impact le plus visible : la modification du paysage à travers la modification de la géométrie des cours d’eau, du relief et des dépôts sédimentaires.
• Enfin ce séisme a tué 28 personnes.

En effet le tremblement de terre a causé un glissement de terrain qui a déplacé 80 millions de tonnes de matière. Le lac a été rempli de boue et de roche, ce qui transforma de manière pérène l'evironnement.

Au niveau économique, les conséquences n'étaient pas si graves. En effet, mis à part la réparation du barrage, et la fermeture du parc pendant plusieurs semaines, ce séisme n'a pas couté "cher" au Wyoming. Le parc a du encaisser certains couts, mais qui ont vite été comblés par l'affluence des Américains dès son ouverture.

Cependant, le tremblement de terre a durement affecté la biodiversité. Ce parc accueille un nombre colossal d'espèces animales différentes, et beaucoup d'animaux ont du mourir à cause de ce séisme. La période de reproduction de la majorité des animaux étant vers le printemps, ce séisme (survenu en aout) a du tuer un grand nombre de bébés ce qui a mis un frein à la reproduction des espèces de ce parc.

Très peu de gens vivent proche de ce parc, il se situe en pleine forêt. Le nombre de mort est assez faible pour un séisme de cette magnitude. On a détecté très peu voir aucun impact sur l'affluence de ce parc à sa réouverture. Bien sûr, certaines maisons proches de lépicentre ont été endommagées, ainsi que certaines routes, etc. Cependant, les dégats restent minimes, ce séisme représente un cout de reconstruction faible.

Ainsi, ce tremblement de terre a transformé le paysage de YellowStone. Il a eu un impact économique et social assez faible, bien que 28 personnes soient décédées. Cependant, l'impact environnemental a été énorme, notamment à travers la transformation des sols, et les conséquences sur la biodiversité de la faune et la flore.

Calcul de la population touchée par chaque niveau d’intensité macrosismique

Nous allons d'abord étudier le jeu de donnée qui contient les coordonnées des lignes d'intensité macrosismique (MMI). Ensuite, nous allons importer la base de données de la densité de population mondiale par la NASA. Cette base de donnée, une fois transformée en coordonées spatiales, sera "superposée" à la carte initiale. On pourra estimer l'apprtenance de points de densité à un niveau d'intensité, en fonction de la coordonnée de ce point, par rapport aux coordonnées des points des lignes d'intensité macrosismique.

I - Imports

#json et os pour analyser le premier dataset avec les intensité macrosismiques
import json
import os

#folium pour créer des cartes
import folium
from folium import plugins
from folium.plugins import HeatMap

#rasterio pour jouer avec la carte de densité
import rasterio as rio
from rasterio.plot import show
from rasterio.warp import calculate_default_transform, reproject, Resampling

#Geopandas pour les cartes
import geopandas

#Pil pour afficher les images de carte
from PIL import Image

#Gdal pour l'étude des coordonnées
from osgeo import gdal, osr

#autres
import numpy as np
import pandas as pd
import math
from math import *
import warnings

from affine import Affine
from pyproj import Proj, transform, CRS #projection de point en coordonnées

import matplotlib.pyplot as plt

II - Analyse du premier set de donnée : courbes d'intensité macrosismiques

fin = pd.DataFrame()

#Ici, on va nettoyer le dataset pour récupérer un dataframe pandas, plus facile à manipuler
#Chaque colonne correspond à une couleur et donc à une courbe d'intensité sismique
#Dans chaque colonne, on a les points qui créent cette courbe
#La 1e colonne correspond à la plus faible intensité, et on monte d'intensité et se déplaçant sur le dataframe

# Opening JSON file 
with open('cont_mmi.json') as json_file: 
    data = json.load(json_file) 
  
    # Print the data of dictionary 
    #'type', 'crs', 'features', 'metadata', 'bbox'])
    for feature in data['features'] :
        df = pd.DataFrame()
        df0 = pd.DataFrame()
        Color = feature['properties']['color']
        Coordinates = feature['geometry']['coordinates']  
        df[Color] = Coordinates
        
        df0 = pd.concat([df0,df], axis=1)

        length = [len(df0[Color][i]) for i in range(len(df[Color]))] #taille des données de chaque intensité
        df_ini=[]
        for i in range(0, len(length)) :
            k = length[i]
            for j in range(int(k)):
                df_ini.append(df0[Color][i][j]) #séparation des coordonées pour avoir un dataframe lisible
        la=len(df_ini)
        df_ini.extend(['']*(1094-la)) #normalisation du dataframe

        fin[Color] = df_ini

fin.head()

	#a0e5ff	#90f2ff	#80ffff	#7cffc7	#7cff90	#bbff4a	#ffff00	#ffe200	#ffc600	#ffaa00	#ff9100
0	[-115.583333, 47.773397]	[-113.971397, 40.966628]	[-111.88948, 46.616863]	[-112.471844, 43.020934]	[-111.090808, 45.60502]	[-111.140725, 43.708319]	[-111.007614, 44.057417]	[-110.758029, 44.253234]	[-110.941058, 44.341171]	[-110.841224, 44.473263]	[-110.874502, 44.537804]
1	[-115.550055, 47.758954]	[-113.955449, 41.066395]	[-111.989313, 46.612945]	[-112.521761, 43.015835]	[-111.140725, 45.590223]	[-111.207281, 43.723548]	[-111.05753, 44.056113]	[-110.90778, 44.276401]	[-110.990975, 44.342224]	[-110.957697, 44.469825]	[-110.924419, 44.538506]
2	[-115.516777, 47.77304]	[-113.897447, 41.216045]	[-111.997739, 46.620085]	[-112.555039, 43.064275]	[-111.209356, 45.622416]	[-111.307115, 43.86213]	[-111.090808, 44.07304]	[-111.007614, 44.263549]	[-111.074169, 44.390799]	[-111.05753, 44.520579]	[-111.007614, 44.602258]
3	[-115.450221, 47.755391]	[-113.896138, 41.249301]	[-111.980709, 46.686597]	[-112.604956, 43.090459]	[-111.291171, 45.738811]	[-111.357032, 43.886585]	[-111.157364, 44.080022]	[-111.074169, 44.270847]	[-111.190642, 44.50224]	[-111.141524, 44.641375]	[-111.015756, 44.641375]
4	[-115.433583, 47.722971]	[-113.965918, 41.415579]	[-111.922757, 46.674481]	[-112.660288, 43.194755]	[-111.317779, 45.805322]	[-111.573339, 43.904352]	[-111.257198, 44.065532]	[-111.207281, 44.333767]	[-111.22392, 44.551418]	[-111.144203, 44.691259]	[-110.990975, 44.683275]

#Visualisons la répartition des points pour quelques niveaux d'intensité
fin1 = fin['#ffff00'][:1093]
fin2 = fin['#80ffff'][:1093]
fin3 = fin['#90f2ff'][:1093]

latitude = [fin1[i][0] for i in range(50)]
longitude = [fin1[i][1] for i in range(50)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude'] #on récupère ici la liste des latitudes et la liste des longitudes

latitude2 = [fin2[i][0] for i in range(525)]
longitude2 = [fin2[i][1] for i in range(525)]
lat_long2 = pd.DataFrame([latitude2, longitude2]).transpose()
lat_long2.columns = ['latitude', 'longitude'] #on récupère ici la liste des latitudes et la liste des longitudes

latitude3 = [fin3[i][0] for i in range(1093)]
longitude3 = [fin3[i][1] for i in range(1093)]
lat_long3 = pd.DataFrame([latitude3, longitude3]).transpose()
lat_long3.columns = ['latitude', 'longitude'] #on récupère ici la liste des latitudes et la liste des longitudes

#On crée une carte pour visualiser ces niveaux d'intensité avec la librairie folium

gdf = geopandas.GeoDataFrame(fin3, geometry=geopandas.points_from_xy(lat_long3.latitude, lat_long3.longitude))
m = folium.Map()

#On adapte l'échelle, avec nos données, en utilisant les latitudes min et max
minx, miny, maxx, maxy = gdf.total_bounds
m.fit_bounds([[miny, minx], [maxy, maxx]])

#On insère des points (petites cercles) à chaque donnée
for i in range(0,1093):
    folium.Circle(
        radius=10,
        location=[lat_long3.longitude[i], lat_long3.latitude[i]],
        popup='The Waterfront',
        color='green',
        fill=False,
    ).add_to(m)

for i in range(0,50):
    folium.Circle(
        radius=10,
        location=[lat_long.longitude[i], lat_long.latitude[i]],
        popup='The Waterfront',
        color='crimson',
        fill=False,
    ).add_to(m)

for i in range(0,525):
    folium.Circle(
        radius=10,
        location=[lat_long2.longitude[i], lat_long2.latitude[i]],
        popup='The Waterfront',
        color='#3186cc',
        fill=False,
    ).add_to(m)
    
m

Make this Notebook Trusted to load map: File -> Trust Notebook

#On va rajouter l'épicentre, qui d'après le premier doc fourni, a les coordonées : 44.7, -111.2

folium.Marker(
    location=[44.7, -111.2],
    popup='Epicentre',
    icon=folium.Icon(color='blue', icon='info-sign')
).add_to(m)

folium.CircleMarker(
    location=[44.7, -111.2],#et on rajoute un cercle sur la carte pour bien le visualiser
    radius=5,
    popup='Laurelhurst Park',
    color='#3186cc',
    fill=True,
    fill_color='#3186cc' #le choix des couleurs est fait par moi même pour accentuer le contraste
).add_to(m)

m

Make this Notebook Trusted to load map: File -> Trust Notebook

#On va ici créer une heatmap de ces mêmes lignes d'intensité macrosismiques,
#dans la continuation de l'étude du dataset, pour visualiser d'une autre façon les données

map_hooray = folium.Map(location=[44.7, -111.2],
                    zoom_start = 7) 

#ici on crée une liste avec la taille des listes de donnée effective pour chaque couleur (colonne)
length_col=[]
for col in fin.columns :
    for i in range(1093) :
        if fin[col][i]== "" :
            length_col.append(i)
            break
        if i==1092 :
            length_col.append(i) 
            
# Liste en compréhension 
j=0
for col in fin.columns :
    heat_data = [[ fin[col][i][1], fin[col][i][0] ] for i in range(length_col[j])]
    HeatMap(heat_data).add_to(map_hooray)
    j+=1

# Display the map
map_hooray

Make this Notebook Trusted to load map: File -> Trust Notebook

III - Analyse de la base de donnée de densité

#On ouvre l'image, et on réucpère sa taille, sa résolution
im = Image.open('pop_density.tif')

width, height = im.size
print("width :",width, "   height :",height)

width : 8640    height : 4320

#On va ici utiliser rasterio pour avoir une analyse statistique des données, notamment la taille des pixels

dataset = gdal.Open('pop_density.tif', gdal.GA_ReadOnly) #ouverture du jeu de données

print("Size is {} x {} x {}".format(dataset.RasterXSize,dataset.RasterYSize,dataset.RasterCount))

geotransform = dataset.GetGeoTransform() #utilisation de la fonction GetGeoTransform() du package
if geotransform:
    print("Origin = ({}, {})".format(geotransform[0], geotransform[3]))
    print("Pixel Size = ({}, {})".format(geotransform[1], geotransform[5]))

Size is 8640 x 4320 x 1
Origin = (-180.0, 89.99999999999994)
Pixel Size = (0.0416666666666667, -0.0416666666666667)

#L'image est sous forme de pixel, ainsi, nous avons besoin de transdformer ces pixels en coordonées que l'on pourra
#ensuite comparer aux coordonées du premier dataset

def pixel_to_coordinates(x, y):
    """Prend en entrée les coordonées d'un pixel (x,y)
       et sort les coordonnées spatiales du pixel correspondant
    """
    
    xoffset, px_w, rot1, yoffset, px_h, rot2 = ds.GetGeoTransform()

    # Les coordonées du pixel sont (x,y)
    # On utilise cette formule pour récupérer les coordonées spatiales de (x,y)
    posX = px_w * x + rot1 * y + xoffset
    posY = rot2 * x + px_h * y + yoffset

    # Décalage au centre du pixel
    posX += px_w / 2.0
    posY += px_h / 2.0


    # Récupère le CRS du dataset, avec la librairire osgeo 
    crs = osr.SpatialReference()
    crs.ImportFromWkt(ds.GetProjectionRef())
    # On créée les latitudes et longitudes CRS avec la ref WGS84
    crsGeo = osr.SpatialReference()
    crsGeo.ImportFromEPSG(4326) # 4326 est l'id EPSG de la latitude/longitude CRS
    t = osr.CoordinateTransformation(crs, crsGeo)
    (long, lat, z) = t.TransformPoint(posX, posY) #ON récupère long et lat, les longitudfes et latitudes
    
    return (lat, long) 

#exemple
pixel_to_coordinates(30, 70)

def is_in_zone(x, y) :
    """Cette fonction permet de définir si un point de coordonée (x,y) appartient à la zone prédéfinit
    Cette zone correspond à la différence entre la surface générée par la courbe d'intensité macrosismique N -
    (grand cercle) moins la surface crée par la courbe d'intensité (N-1) - (petit cercel)
    """
    if (x < miny_g or x > maxy_g or y < minx_g or y > maxx_g) :
        return False #vérifie si le point est à l'intérieur du "grand cercle" 
    
    if (x<maxy_p and x >miny_p) :
        if (y<maxx_p and y >minx_p) :
            return False #vérifie si le point est à l'exterieur du "petit cercle" 
    
    return True

#exemple
is_in_zone(42.49999999999991, -115.5)

False

img = Image.open('pop_density.tif')
img.close()

fname = 'pop_density.tif'
warnings.filterwarnings("ignore", category=DeprecationWarning)

# Analyse l'image avec rasterio
with rasterio.open(fname, 'r+') as r:
    T0 = r.transform #Transformation afffine du pixel du coin en haut à gauche de l'image
    p1 = Proj(r.crs) #Projection avec crs
    A = r.read()  #Valeurs des pixels

#Créée toutes les colonnes et lignes
cols, rows = np.meshgrid(np.arange(A.shape[2]), np.arange(A.shape[1]))

#Récupère la transformation affine pour les centres des pixels
T1 = T0 * Affine.translation(0.5, 0.5)
# Convertit les indices des colonnes/lignes des pixels en abscisses et ordonnées
rc2en = lambda r, c: T1 * (c, r)

# All eastings and northings (abscisses et ordonnées)
eastings, northings = np.vectorize(rc2en, otypes=[np.float, np.float])(rows, cols)

# Projection de toutes les latitudes et longitudes
p2 = Proj(proj='latlong',datum='WGS84')
lats, longs = transform(p1, p2, eastings, northings) #récupère toutes les latitudes/longitudes

def nombre_de_points(petit_cercle, grand_cercle) :
    """Récupère, le nombre de points de densité dans la surface entre le "grand" cercle et le "petit" cercle
    """
    pixel = im.load()
    Nombre_prs=0

    for x in range(height):
        for y in range(width): #on parcours tous les pixels
            lat, long = lats[x, y], longs[x, y] #On récupère les latitude et longitudes
            if is_in_zone(lat, long) : #On ajoute au compteur si les coordonnées sont dans la zone donnée
                Nombre_prs+=1 

    return int(Nombre_prs)

IV - Analyse par niveau d'intensité

Ici, pour chaque niveau d'intensité, on va récupérer la superficie de la surface qui est analysée, la densité de population pour cette surface, et enfin le nombre d'habitants estimés en fonction des ces deux dernières données

#Niveau 1 - De plus forte intensité macrosismique

#Ici, puisque c'est le premier cercle, il n'y a pas de cercle plus petit, donc on met les coordonées de l'épicentre
minx_p, miny_p, maxx_p, maxy_p = -111.2, 44.70, -111.23, 44.71 #épicentre

#A chaque fois on définit le grand cercle, on récupère les données du grand cercle, et ses coordonnées limites
grand_cercle = fin['#ff9100'][:1093] #cercle petit 

latitude = [grand_cercle[i][0] for i in range(11)]
longitude = [grand_cercle[i][1] for i in range(11)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_grand = geopandas.GeoDataFrame(grand_cercle[:11], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_grand.total_bounds


#on déinit une superficie en fonction des coordonées limites
Superficie_1 = abs(minx_g-maxx_g)*abs(maxy_g-miny_g)*5000 
print(int(Superficie_1), "mètres carrés, est le superficie de la zone associée en mètres")

#On définit la densité de cette zone comme le nombre de points de densité qui appartiennent à cette zone
Densite_1 = nombre_de_points(0, grand_cercle)
print("On récupère", Densite_1, "points de densité")

Habitants_1 = int(Superficie_1)* Densite_1 / 100000
print("Sur cette zone, il y a donc", floor(Habitants_1), "habitants")

184 mètres carrés, est le superficie de la zone associée en mètres
On récupère 20 points de densité
Sur cette zone, il y a donc 0 habitants

#Niveau 2 - celui d'intensité sismique qui suit

#A nouveau on définit le petit cercle, on récupère les données du grand cercle, et ses coordonnées limites
petit_cercle = fin['#ff9100'][:1093]#cercle petit  

latitude = [petit_cercle[i][0] for i in range(11)]
longitude = [petit_cercle[i][1] for i in range(11)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_petit = geopandas.GeoDataFrame(petit_cercle[:11], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_p, miny_p, maxx_p, maxy_p = gdf_petit.total_bounds #on récupère les coordonées limites


#A nouveau on définit le petit cercle, on récupère les données du grand cercle, et ses coordonnées limites
grand_cercle = fin['#ffaa00'][:1093] #cercle petit 

latitude = [grand_cercle[i][0] for i in range(14)]
longitude = [grand_cercle[i][1] for i in range(14)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_grand = geopandas.GeoDataFrame(grand_cercle[:14], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_grand.total_bounds#on récupère les coordonées limites


#Puis calcul de la superficie, densité et population
Superficie_2 = abs(minx_g-maxx_g)*abs(maxy_g-miny_g)*5000 - Superficie_1 #superficie de la zone
print(int(Superficie_2), "mètres carrés, est le superficie de la zone associée en mètres")

Densite_2 = nombre_de_points(petit_cercle, grand_cercle)
print("On récupère", Densite_2, "points de densité")

Habitants_2 = int(Superficie_2)* Densite_2 / 100000
print("Sur cette zone, il y a donc", floor(Habitants_2), "habitants")

565 mètres carrés, est le superficie de la zone associée en mètres
On récupère 68 points de densité
Sur cette zone, il y a donc 0 habitants

#Niveau 3

#A nouveau on définit le petit cercle, on récupère les données du grand cercle, et ses coordonnées limites
petit_cercle = fin['#ffaa00'][:1093] #cercle petit 

latitude = [petit_cercle[i][0] for i in range(14)]
longitude = [petit_cercle[i][1] for i in range(14)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_petit = geopandas.GeoDataFrame(petit_cercle[:14], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_petit.total_bounds

#A nouveau on définit le grand cercle, on récupère les données du grand cercle, et ses coordonnées limites
grand_cercle = fin['#ffc600'][:1093] #cercle petit 

latitude = [grand_cercle[i][0] for i in range(19)]
longitude = [grand_cercle[i][1] for i in range(19)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_grand = geopandas.GeoDataFrame(grand_cercle[:19], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_grand.total_bounds


#Puis calcul de la superficie, densité et population
Superficie_3 = abs(minx_g-maxx_g)*abs(maxy_g-miny_g)*5000 - Superficie_2 #superficie de la zone
print(int(Superficie_3), "mètres carrés, est le superficie de la zone associée en mètres")

Densite_3 = nombre_de_points(petit_cercle, grand_cercle)
print("On récupère", Densite_3, "points de densité")

Habitants_3 = int(Superficie_3)* Densite_3 / 100000
print("Sur cette zone, il y a donc", floor(Habitants_3), "habitants")

1221 mètres carrés, est le superficie de la zone associée en mètres
On récupère 184 points de densité
Sur cette zone, il y a donc 2 habitants

#Niveau 4

#A nouveau on définit le petit cercle, on récupère les données du grand cercle, et ses coordonnées limites
petit_cercle = fin['#ffc600'][:1093] #cercle petit 

latitude = [petit_cercle[i][0] for i in range(19)]
longitude = [petit_cercle[i][1] for i in range(19)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_petit = geopandas.GeoDataFrame(petit_cercle[:19], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_petit.total_bounds

#A nouveau on définit le grand cercle, on récupère les données du grand cercle, et ses coordonnées limites
grand_cercle = fin['#ffe200'][:1093] #cercle petit 

latitude = [grand_cercle[i][0] for i in range(28)]
longitude = [grand_cercle[i][1] for i in range(28)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_grand = geopandas.GeoDataFrame(grand_cercle[:28], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_grand.total_bounds


#Puis calcul de la superficie, densité et population
Superficie_4 = abs(minx_g-maxx_g)*abs(maxy_g-miny_g)*5000 - Superficie_3 #superficie de la zone
print(int(Superficie_4), "mètres carrés, est le superficie de la zone associée en mètres")

Densite_4 = nombre_de_points(petit_cercle, grand_cercle)
print("On récupère", Densite_4, "points de densité")

Habitants_4 = int(Superficie_4)* Densite_4 / 100000
print("Sur cette zone, il y a donc", floor(Habitants_4), "habitants")

2339 mètres carrés, est le superficie de la zone associée en mètres
On récupère 405 points de densité
Sur cette zone, il y a donc 9 habitants

#Niveau 5

#A nouveau on définit le petit cercle, on récupère les données du grand cercle, et ses coordonnées limites
petit_cercle = fin['#ffe200'][:1093] #cercle petit 

latitude = [petit_cercle[i][0] for i in range(28)]
longitude = [petit_cercle[i][1] for i in range(28)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_petit = geopandas.GeoDataFrame(petit_cercle[:28], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_petit.total_bounds

#A nouveau on définit le grand cercle, on récupère les données du grand cercle, et ses coordonnées limites
grand_cercle = fin['#ffff00'][:1093] #cercle petit 

latitude = [grand_cercle[i][0] for i in range(50)]
longitude = [grand_cercle[i][1] for i in range(50)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_grand = geopandas.GeoDataFrame(grand_cercle[:50], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_grand.total_bounds


#Puis calcul de la superficie, densité et population
Superficie_5 = abs(minx_g-maxx_g)*abs(maxy_g-miny_g)*5000 - Superficie_4 #superficie de la zone
print(int(Superficie_5), "mètres carrés, est le superficie de la zone associée en mètres")

Densite_5 = nombre_de_points(petit_cercle, grand_cercle)
print("On récupère", Densite_5, "points de densité")

Habitants_5 = int(Superficie_5)* Densite_5 / 100000
print("Sur cette zone, il y a donc", floor(Habitants_5), "habitants")

5401 mètres carrés, est le superficie de la zone associée en mètres
On récupère 892 points de densité
Sur cette zone, il y a donc 48 habitants

#Niveau 6

#A nouveau on définit le petit cercle, on récupère les données du grand cercle, et ses coordonnées limites
petit_cercle = fin['#ffff00'][:1093] #cercle petit 

latitude = [petit_cercle[i][0] for i in range(50)]
longitude = [petit_cercle[i][1] for i in range(50)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_petit = geopandas.GeoDataFrame(petit_cercle[:50], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_petit.total_bounds

#A nouveau on définit le grand cercle, on récupère les données du grand cercle, et ses coordonnées limites
grand_cercle = fin['#bbff4a'][:1093] #cercle petit 

latitude = [grand_cercle[i][0] for i in range(67)]
longitude = [grand_cercle[i][1] for i in range(67)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_grand = geopandas.GeoDataFrame(grand_cercle[:67], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_grand.total_bounds


#Puis calcul de la superficie, densité et population
Superficie_6 = abs(minx_g-maxx_g)*abs(maxy_g-miny_g)*5000 - Superficie_5 #superficie de la zone
print(int(Superficie_6), "mètres carrés, est le superficie de la zone associée en mètres")

Densite_6 = nombre_de_points(petit_cercle, grand_cercle)
print("On récupère", Densite_6, "points de densité")

Habitants_6 = int(Superficie_6)* Densite_6 / 100000
print("Sur cette zone, il y a donc", floor(Habitants_6), "habitants")

9034 mètres carrés, est le superficie de la zone associée en mètres
On récupère 1636 points de densité
Sur cette zone, il y a donc 147 habitants

#Niveau 7

#A nouveau on définit le petit cercle, on récupère les données du grand cercle, et ses coordonnées limites
petit_cercle = fin['#bbff4a'][:1093] #cercle petit 

latitude = [petit_cercle[i][0] for i in range(67)]
longitude = [petit_cercle[i][1] for i in range(67)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_petit = geopandas.GeoDataFrame(petit_cercle[:67], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_petit.total_bounds

#A nouveau on définit le grand cercle, on récupère les données du grand cercle, et ses coordonnées limites
grand_cercle = fin['#7cff90'][:1093] #cercle petit 

latitude = [grand_cercle[i][0] for i in range(130)]
longitude = [grand_cercle[i][1] for i in range(130)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_grand = geopandas.GeoDataFrame(grand_cercle[:130], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_grand.total_bounds


#Puis calcul de la superficie, densité et population
Superficie_7 = abs(minx_g-maxx_g)*abs(maxy_g-miny_g)*5000 - Superficie_6 #superficie de la zone
print(int(Superficie_7), "mètres carrés, est le superficie de la zone associée en mètres")

Densite_7 = nombre_de_points(petit_cercle, grand_cercle)
print("On récupère", Densite_7, "points de densité")

Habitants_7 = int(Superficie_7)* Densite_7 / 100000
print("Sur cette zone, il y a donc", floor(Habitants_7), "habitants")

27923 mètres carrés, est le superficie de la zone associée en mètres
On récupère 4240 points de densité
Sur cette zone, il y a donc 1183 habitants

#Niveau 8

#A nouveau on définit le petit cercle, on récupère les données du grand cercle, et ses coordonnées limites
petit_cercle = fin['#7cff90'][:1093] #cercle petit 

latitude = [petit_cercle[i][0] for i in range(130)]
longitude = [petit_cercle[i][1] for i in range(130)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_petit = geopandas.GeoDataFrame(petit_cercle[:130], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_petit.total_bounds

#A nouveau on définit le grand cercle, on récupère les données du grand cercle, et ses coordonnées limites
grand_cercle = fin['#7cffc7'][:1093] #cercle petit 

latitude = [grand_cercle[i][0] for i in range(257)]
longitude = [grand_cercle[i][1] for i in range(257)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_grand = geopandas.GeoDataFrame(grand_cercle[:257], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_grand.total_bounds


#Puis calcul de la superficie, densité et population
Superficie_8 = abs(minx_g-maxx_g)*abs(maxy_g-miny_g)*5000 - Superficie_7 #superficie de la zone
print(int(Superficie_8), "mètres carrés, est le superficie de la zone associée en mètres")

Densite_8 = nombre_de_points(petit_cercle, grand_cercle)
print("On récupère", Densite_8, "points de densité")

Habitants_8 = int(Superficie_8)* Densite_8 / 100000
print("Sur cette zone, il y a donc", floor(Habitants_8), "habitants")

36348 mètres carrés, est le superficie de la zone associée en mètres
On récupère 7330 points de densité
Sur cette zone, il y a donc 2664 habitants

#Niveau 9

#A nouveau on définit le petit cercle, on récupère les données du grand cercle, et ses coordonnées limites
petit_cercle = fin['#7cffc7'][:1093] #cercle petit 

latitude = [petit_cercle[i][0] for i in range(257)]
longitude = [petit_cercle[i][1] for i in range(257)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_petit = geopandas.GeoDataFrame(petit_cercle[:257], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_petit.total_bounds

#A nouveau on définit le grand cercle, on récupère les données du grand cercle, et ses coordonnées limites
grand_cercle = fin['#80ffff'][:1093] #cercle petit 

latitude = [grand_cercle[i][0] for i in range(526)]
longitude = [grand_cercle[i][1] for i in range(526)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_grand = geopandas.GeoDataFrame(grand_cercle[:526], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_grand.total_bounds


#Puis calcul de la superficie, densité et population
Superficie_9 = abs(minx_g-maxx_g)*abs(maxy_g-miny_g)*5000 - Superficie_8 #superficie de la zone
print(int(Superficie_9), "mètres carrés, est le superficie de la zone associée en mètres")

Densite_9 = nombre_de_points(petit_cercle, grand_cercle)
print("On récupère", Densite_9, "points de densité")

Habitants_9 = int(Superficie_9)* Densite_9 / 100000
print("Sur cette zone, il y a donc", floor(Habitants_9), "habitants")

109141 mètres carrés, est le superficie de la zone associée en mètres
On récupère 16624 points de densité
Sur cette zone, il y a donc 18143 habitants

#Niveau 10

#A nouveau on définit le petit cercle, on récupère les données du grand cercle, et ses coordonnées limites
petit_cercle = fin['#80ffff'][:1093] #cercle petit 

latitude = [petit_cercle[i][0] for i in range(526)]
longitude = [petit_cercle[i][1] for i in range(526)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_petit = geopandas.GeoDataFrame(petit_cercle[:526], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_petit.total_bounds

#A nouveau on définit le grand cercle, on récupère les données du grand cercle, et ses coordonnées limites
grand_cercle = fin['#90f2ff'][:1093] #cercle petit 

latitude = [grand_cercle[i][0] for i in range(1093)]
longitude = [grand_cercle[i][1] for i in range(1093)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_grand = geopandas.GeoDataFrame(grand_cercle[:1093], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_grand.total_bounds


#Puis calcul de la superficie, densité et population
Superficie_10 = abs(minx_g-maxx_g)*abs(maxy_g-miny_g)*5000 - Superficie_9 #superficie de la zone
print(int(Superficie_10), "mètres carrés, est le superficie de la zone associée en mètres")

Densite_10 = nombre_de_points(petit_cercle, grand_cercle)
print("On récupère", Densite_10, "points de densité")

Habitants_10 = int(Superficie_10)* Densite_10 / 100000
print("Sur cette zone, il y a donc", floor(Habitants_10), "habitants")

221890 mètres carrés, est le superficie de la zone associée en mètres
On récupère 38113 points de densité
Sur cette zone, il y a donc 84568 habitants

#Niveau 11

#A nouveau on définit le petit cercle, on récupère les données du grand cercle, et ses coordonnées limites
petit_cercle = fin['#90f2ff'][:1093] #cercle petit 

latitude = [petit_cercle[i][0] for i in range(1093)]
longitude = [petit_cercle[i][1] for i in range(1093)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_petit = geopandas.GeoDataFrame(petit_cercle[:1093], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_petit.total_bounds

#A nouveau on définit le grand cercle, on récupère les données du grand cercle, et ses coordonnées limites
grand_cercle = fin['#a0e5ff'][:1093] #cercle petit

latitude = [grand_cercle[i][0] for i in range(140)]
longitude = [grand_cercle[i][1] for i in range(140)]
lat_long = pd.DataFrame([latitude, longitude]).transpose()
lat_long.columns = ['latitude', 'longitude']

gdf_grand = geopandas.GeoDataFrame(grand_cercle[:140], geometry=geopandas.points_from_xy(lat_long.latitude, lat_long.longitude))
minx_g, miny_g, maxx_g, maxy_g = gdf_grand.total_bounds


#Puis calcul de la superficie, densité et population
Superficie_11 = abs(minx_g-maxx_g)*abs(maxy_g-miny_g)*5000 - Superficie_10 #superficie de la zone
print(int(Superficie_11), "mètres carrés, est le superficie de la zone associée en mètres")

Densite_11 = nombre_de_points(petit_cercle, grand_cercle)
print("On récupère", Densite_11, "points de densité")

Habitants_11 = int(Superficie_11)* Densite_11 / 100000
print("Sur cette zone, il y a donc", floor(Habitants_11), "habitants")

133945 mètres carrés, est le superficie de la zone associée en mètres
On récupère 41020 points de densité
Sur cette zone, il y a donc 54944 habitants

#BILAN

Habitants_tot = Habitants_1 + Habitants_2 + Habitants_3 + Habitants_4 + Habitants_5 + Habitants_6 + Habitants_7 + Habitants_8 + Habitants_9 + Habitants_10 + Habitants_11
print("Il y a à peu près", floor(Habitants_tot), "habitants sur l'ensemble des zones touchées par le séisme")

print("")

Superficie_tot = Superficie_1 + Superficie_2 + Superficie_3 + Superficie_4 + Superficie_5 + Superficie_6 + Superficie_7 + Superficie_8 + Superficie_9 + Superficie_10 + Superficie_11
print("Une surface totale de", floor(Superficie_tot), "mètres carrés a été affectée par le séisme, soit", round(Superficie_tot*1e-6, 2), "en kilomètres carrés")

Il y a à peu près 161713 habitants sur l'ensemble des zones touchées par le séisme

Une surface totale de 547996 mètres carrés a été affectée par le séisme, soit 0.55 en kilomètres carrés